{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 使用K-Means算法进行聚类"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## K-Means算法基本思想\n",
    "\n",
    "事先确定常数K，常数K意味着最终的聚类类别数，首先随机选定初始点为质心，并通过计算每一个样本与质心之间的相似度(这里为欧式距离)，将样本点归到最相似的类中，接着，重新计算每个类的质心(即为类中心)，重复这样的过程，直到质心不再改变，最终就确定了每个样本所属的类别以及每个类的质心。由于每次都要计算所有的样本与每一个质心之间的相似度，故在大规模的数据集上，K-Means算法的收敛速度比较慢。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## K-Means算法步骤\n",
    "\n",
    "（1）任选K个初始聚类中心：$ Z_{1}(1), \\quad Z_{2}(1), \\quad \\ldots, \\quad Z_{K}(1) $。括号内序号：迭代运算的次序号。\n",
    "\n",
    "（2）按最小距离原则将其余样本分配到K个聚类中心中的某一个，即：若\n",
    "$ \\min \\left\\{\\left\\|\\boldsymbol{X}-\\boldsymbol{Z}_{i}(k)\\right\\|, i=1,2, \\cdots, K\\right\\}=\\left\\|\\boldsymbol{X}-\\boldsymbol{Z}_{j}(k)\\right\\|=D_{j}(k) $，则$ X \\in S_{j}(k) $。其中k为迭代运算次序号，K为聚类中心的个数。\n",
    "\n",
    "（3）计算各个聚类中心的新向量值：$ Z_{j}(k+1) \\quad j=1,2, \\cdots, K $\n",
    "$ \\boldsymbol{Z}_{j}(k+1)=\\frac{1}{N_{j}} \\sum_{X \\in S_{j}(k)} \\boldsymbol{X} \\quad j=1,2, \\cdots, K , \\boldsymbol{N}_{j} :$第j类的样本数\n",
    "\n",
    "（4）如果$ \\boldsymbol{Z}_{j}(k+1) \\neq \\boldsymbol{Z}_{j}(k) \\quad j=1,2, \\cdots, K $，则回到（2），将模式样本逐个重新分类，重复迭代计算。\n",
    "如果$\\boldsymbol{Z}_{j}(k+1) = \\boldsymbol{Z}_{j}(k) \\quad j=1,2, \\cdots, K$，算法收敛，计算完毕。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "使用K-Means算法实现对样本的分类，选择的样本为二维模式，设置两个聚类中心，画出样本聚类情况，并判断[2,3],[6,9]分别属于哪一类。样本如下：\n",
    " \n",
    "        [0,0], [1,0], [0,1], [1,1], [2,1], [1,2], [2,2], [3,2], [6,6], [7,6],\n",
    "        [8,6], [6,7], [7,7], [8,7], [9,7], [7,8], [8,8], [9,8], [8,9], [9,9]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "第1次迭代\n",
      "中心点: {0: array([0, 0]), 1: array([1, 0])}\n",
      "分组情况:\n",
      "{   0: [array([0, 0]), array([0, 1])],\n",
      "    1: [   array([1, 0]), array([1, 1]), array([2, 1]), array([1, 2]), array([2, 2]), array([3, 2]), array([6, 6]),\n",
      "           array([7, 6]), array([8, 6]), array([6, 7]), array([7, 7]), array([8, 7]), array([9, 7]), array([7, 8]),\n",
      "           array([8, 8]), array([9, 8]), array([8, 9]), array([9, 9])]}\n",
      "----------------------------------------------------------------------\n",
      "第2次迭代\n",
      "中心点: {0: array([0. , 0.5]), 1: array([5.66666667, 5.33333333])}\n",
      "分组情况:\n",
      "{   0: [   array([0, 0]), array([1, 0]), array([0, 1]), array([1, 1]), array([2, 1]), array([1, 2]), array([2, 2]),\n",
      "           array([3, 2])],\n",
      "    1: [   array([6, 6]), array([7, 6]), array([8, 6]), array([6, 7]), array([7, 7]), array([8, 7]), array([9, 7]),\n",
      "           array([7, 8]), array([8, 8]), array([9, 8]), array([8, 9]), array([9, 9])]}\n",
      "----------------------------------------------------------------------\n",
      "第3次迭代\n",
      "中心点: {0: array([1.25 , 1.125]), 1: array([7.66666667, 7.33333333])}\n",
      "分组情况:\n",
      "{   0: [   array([0, 0]), array([1, 0]), array([0, 1]), array([1, 1]), array([2, 1]), array([1, 2]), array([2, 2]),\n",
      "           array([3, 2])],\n",
      "    1: [   array([6, 6]), array([7, 6]), array([8, 6]), array([6, 7]), array([7, 7]), array([8, 7]), array([9, 7]),\n",
      "           array([7, 8]), array([8, 8]), array([9, 8]), array([8, 9]), array([9, 9])]}\n",
      "----------------------------------------------------------------------\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAD4CAYAAADFAawfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAQg0lEQVR4nO3dfWxd9X3H8c/HcUJwguqqcSsKxA4b6kRpGcRDPJR2Kp3UNX2Ypk3LZKYOqfNWbcWkTFUpfzSblO2PoYZMa6tZdN0fsdpuFA3EGINBS7uxsToElZKwiofYCQ+tQTiscYId+7s/jtPYnh+uw733fO+975dkXZ9zr28+HOGPjn/n3vt1RAgAkFdb2QEAAMujqAEgOYoaAJKjqAEgOYoaAJJrr8WTbtq0KXp6emrx1ADQlPbt2/dKRHQtdl9Nirqnp0fDw8O1eGoAaEq2R5a6j6UPAEiOogaA5ChqAEiOogaA5ChqoMEs/HgePq6n+VHUQAPZuVPaseN0OUcU2zt3lpmqPENDUk+P1NZW3A4NNWcOihpoEBHS+Li0Z8/pst6xo9geH2+9M+uhIam/XxoZKf7bR0aK7XqXdT1yuBYfc9rb2xu8jhqovrnlfMrAgLR7t2SXl6sMPT1FKS7U3S0dOtR4OWzvi4jeRe+jqIHGElH8iX3KzEzrlbRUHIPF6ssujkmj5ViuqFn6ABrIqTPqueauWbeSzZtXt7+Rc1DUQIOYu+wxMFCcrQ0MzF+zbiW7dkkdHfP3dXQU+5stR00+6wNA9dlSZ+f8Nendu4v7Ojtbb/mjr6+4vfVWaXS0OIPdtev0/mbKwRo10GAi5pfywm00JtaogSaysJQp6eZHUQNAchQ1ACRHUQNAchQ1ACRHUQNAchQ1ACRHUQNAchQ1ACRHUQNAchQ1ACRHUQNAchQ1ACRHUQNAchUVte0dtp+y/SPb37C9vtbBAOSWYQJ4hgz1yLFiUds+T9KNknoj4hJJayRtr24MAI0kwwTwDBnqlWPFwQGzRf1fki6V9Lqkf5L01xHxwFI/w+AAoLllmACeIUM1c7ypwQER8YKk2ySNSnpJ0tHFStp2v+1h28NjY2OVpwPQcEZHV7e/WTPUK0clSx9vlfQJSVskvVPSBtvXL3xcRAxGRG9E9HZ1dVUvIYB0MkwAz5ChXjkquZj4IUnPR8RYRExJukvS1dWLAKDRZJgAniFDvXJUUtSjkq603WHbkq6TdLB6EQA0mr4+aXCwWIe1i9vBwfpOAM+QoV45KppCbvvPJP2OpJOS9kv6VES8sdTjuZgIAKuz3MXE9kqeICK+KOmLVU0FAKgI70wEgOQoagBIjqIGkNPJybITpEFRA8jn2CvSbRcVt6CoAST09L3SiXHp6X8uO0kKFDWAfPbPfqLREyV9HF4yFDWAXI6/Jr30RPH9i/ul4+NlpkmhotdRA0DN3HOjdOBuSbNvvps5KbW1S9OTUttaaffFxbYkydK7f0P62J6SwpaDM2oA5br6RmlDlzR1XDpxVJo8Jk1NFPdNHSu2Txwt7t/QJV31mXLzloCiBlCuTb8offo/pMs/Ka09e/HHtJ8tbf196dOPFo9vMRQ1gPK1nyVtu01632eldRvm37d2g3TtZ6WP/JXUvq6cfCWjqAHkcfixYqlDkjy7Lj11TDr83+VlSoCiBpDD5IT0/PeL79vXSxd/rLiVpOe/V6xRtyiKGmgwTTt5+9mHpOk3pI3vkG64T/rtvy9uN76j2P/MQ7XPcIZqniMiqv61devWAFB9e/dGdHREFPOui6+OjmJ/w+d49CsR//DJiBOvz99//Gix/9Gv1D7DGahWDknDsUSnVjQ4YLUYHADURrNN3m70DNXMsdzgAIoaaCBtbcU520K2NDPTWjkyZKhmjuWKmjVqoIG00uTtRshQrxwUNdBAWmnydiNkqFuOpRav38wXFxOB2tm7N6K7O8Iubut98SxTjgwZqpVDXEwEgNxYowaABkZRA0ByFDUAJEdRA0ByFDUAJEdRA0ByFDUAJEdRA0ByFDUAJEdRA0ByFDUAJEdRA0ByFDUAJEdRA0ByFRW17U7bd9p+2vZB21fVOhiAxbXM5O0GyVAP7RU+bo+k+yPit2yvk9Sx0g8AqL6hIam/X5qYKLZHRoptSerra60cGTLUy4qDA2y/RdITki6MCqcMMDgAqI1mm7zd6Bmq6c0ODtgiaUzS123vt32H7Q2L/CP9todtD4+Njb3JyAAWMzq6uv3NnCNDhnqppKjbJV0u6asRcZmkY5I+v/BBETEYEb0R0dvV1VXlmACk1pq83QgZ6qWSoj4i6UhEPDa7faeK4gZQZy01ebsBMtTLikUdES9LOmz7XbO7rpN0oKapACyqr08aHCzWYe3idnCw/hfPMuTIkKFeKppCbvuXJd0haZ2k5yTdEBGvLfV4LiYCwOosdzGxopfnRcQTkhZ9AgBAbfHORABIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIjqIGgOQoagBIruKitr3G9n7b99YyEABgvtWcUQ9IOlirIACAxVVU1LbPl7RN0h21jQMAWKjSM+rbJX1O0sxSD7Ddb3vY9vDY2Fg1sgEAVEFR2/6opJ9GxL7lHhcRgxHRGxG9XV1dVQsIAK2ukjPqayR93PYhSd+U9EHbe2uaCgDwcysWdUTcEhHnR0SPpO2SHo6I62ueDAAgiddRA0B67at5cER8V9J3a5IEALAozqgBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gzilh+G0BLoaiz2blT2rHjdDlHFNs7d5aZCkCJKOpMIqTxcWnPntNlvWNHsT0+zpk10KJW9XnUqDFb2r27+H7PnuJLkgYGiv12edkAlMZRg7O03t7eGB4ervrztowIqW3OHzszM5Q00ORs74uI3sXuY+kjm1PLHXPNXbMG0HIo6kzmrkkPDBRn0gMD89esAbQc1qgzsaXOzvlr0qfWrDs7Wf4AWhRr1BlFzC/lhdsAmg5r1I1mYSlT0kBLo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBILkVi9r2Bba/Y/uA7adsD9QjWGmGhqSenmIUVk9Psd2KGTLlAFpcJYMDTkq6OSIet32OpH22H4yIAzXOVn9DQ1J/vzQxUWyPjBTbktTX1zoZMuUAsPrBAbbvlvQ3EfHgUo9p2MEBPT1FIS3U3S0dOtQ6GTLlAFrEcoMDVlXUtnskfU/SJRHx+oL7+iX1S9LmzZu3jiz2S55dW9vicwntYn5hq2TIlANoEVWZ8GJ7o6RvS7ppYUlLUkQMRkRvRPR2dXWdedoybd68uv3NmiFTDgCVFbXttSpKeigi7qptpBLt2iV1dMzf19FR7G+lDJlyAKjoVR+W9DVJByPiS7WPVKK+PmlwsFiHtYvbwcH6XjzLkCFTDgArr1Hbfp+k70t6UtKpxckvRMR9S/1Mw15MBICSLLdGveLL8yLi3yUxBhsASsI7EwEgOYoaAJKjqAEgOYoaAJKjqAEgOYoaAJKjqM/Q5Ek+7wJAfVDUZ+DVn72hX9n1b3r1Z2+UHQVAC6Coz8ADB36io8en9OCBn5QdBUALoKjPwD8OHy5u9x0pOQmAVkBRr9LRiSk9+cJRSdIPj4zr6PGpkhMBaHaVjOJqabfc9UPd9+TLOvXhVdMzobVr2jQ1Pa21a9p09V8+pDVtxUeh2Na295yrv/jN95QZGUCT4Yx6BX9w7YV624Z1OjE1rddPnNSxyWlNTE5LkiYmp3Vssth/Ympab9u4Tp+6dkvJiQE0mzxFnWXi9YIcFz5wt/7lpmu1/YrNWr928cO1fm2btl+xWfcPvF8Xdm2segamfwOtLcfSR5aJ10vkOEvSn/f1adPGs/S3jzyrY7Nn1JLUsW6N/ugDv6Abr7uophkk8aH9QIvKcUZ9662ni+mUiYlif6Ic+0Ze+3lJr5k9chOT09o/+lrdMgBoPTmKenR0dftLyHF8clr/+dyrkoqljg+/+1ytby8O36PPvqoTU9OL/2wVMwBoTTmKOsvE62VyPPLjMU2enNHbzzlL3+q/Sl/uu1zf+sOr9PZzztIbJ2f0yI/Hap4BQGvKUdRZJl4vk+PF8ePa9t5z9fCf/qouvaBTknTpBZ166OYPaNt7z9ULrx2veQYALSoiqv61devWWLW9eyO6uyPs4nbv3tU/RzVkyJEhA4C6kjQcS3TqilPIzwRTyAFgdZabQp5j6QMAsCSKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSo6gBIDmKGgCSq6iobX/Y9v/Yfsb252uSJMvk7Sw5AGDWilPIba+R9GVJvybpiKQf2L4nIg5ULUWWydtZcgDAHJWcUV8h6ZmIeC4iJiV9U9Inqpoiy+TtLDkAYI5Kivo8SYfnbB+Z3TeP7X7bw7aHx8ZWOeg1y+TtLDkAYI6qXUyMiMGI6I2I3q6urtX9cJbJ21lyAMAclRT1C5IumLN9/uy+6skyeTtLDgCYo5Ki/oGki2xvsb1O0nZJ91Q1RV+fNDgodXdLdnE7OFj/C3hZcgDAHBVNIbf9EUm3S1oj6e8iYtlTTKaQA8DqLDeFfMWX50lSRNwn6b6qpgIAVIR3JgJAchQ1ACRHUQNAchQ1ACRX0as+Vv2k9pikkTP88U2SXqlinEbGsZiP4zEfx+O0ZjgW3RGx6LsFa1LUb4bt4aVeotJqOBbzcTzm43ic1uzHgqUPAEiOogaA5DIW9WDZARLhWMzH8ZiP43FaUx+LdGvUAID5Mp5RAwDmoKgBILk0RV2XAboNwvYFtr9j+4Dtp2wPlJ2pbLbX2N5v+96ys5TNdqftO20/bfug7avKzlQm2ztmf09+ZPsbtteXnanaUhT1nAG6vy7pYkm/a/viclOV6qSkmyPiYklXSvrjFj8ekjQg6WDZIZLYI+n+iPglSZeqhY+L7fMk3SipNyIuUfFRzNvLTVV9KYpa9Rig20Ai4qWIeHz2+/9V8Yv4/+ZUtgrb50vaJumOsrOUzfZbJL1f0tckKSImI2K81FDla5d0tu12SR2SXiw5T9VlKeqKBui2Its9ki6T9FjJUcp0u6TPSZopOUcGWySNSfr67FLQHbY3lB2qLBHxgqTbJI1KeknS0Yh4oNxU1ZelqLEI2xslfVvSTRHxetl5ymD7o5J+GhH7ys6SRLukyyV9NSIuk3RMUste07H9VhV/fW+R9E5JG2xfX26q6stS1LUfoNtgbK9VUdJDEXFX2XlKdI2kj9s+pGJJ7IO295YbqVRHJB2JiFN/Yd2porhb1YckPR8RYxExJekuSVeXnKnqshR17QfoNhDbVrEGeTAivlR2njJFxC0RcX5E9Kj4/+LhiGi6M6ZKRcTLkg7bftfsruskHSgxUtlGJV1pu2P29+Y6NeHF1YpmJtZaRJy0/SeS/lWnB+g+VXKsMl0j6fckPWn7idl9X5idXQl8RtLQ7EnNc5JuKDlPaSLiMdt3Snpcxaul9qsJ307OW8gBILksSx8AgCVQ1ACQHEUNAMlR1ACQHEUNAMlR1ACQHEUNAMn9H+LYD3DXQFfhAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot\n",
    "from pprint import pprint\n",
    "\n",
    "class K_Means(object):\n",
    "    def __init__(self, k=2, tolerance=0.0001, max_iter=300):\n",
    "        \"\"\"\n",
    "        :param k:分组数\n",
    "        :param tolerance:中心点误差\n",
    "        :param max_iter:迭代次数\n",
    "        \"\"\"\n",
    "        self.k_ = k\n",
    "        self.tolerance_ = tolerance\n",
    "        self.max_iter_ = max_iter\n",
    "\n",
    "    def fit(self, data):\n",
    "        \"\"\"  k均值计算 \"\"\"\n",
    "        self.centers_ = {}  # 中心点\n",
    "        for i in range(self.k_):\n",
    "            self.centers_[i] = data[i]\n",
    "\n",
    "        for i in range(self.max_iter_):\n",
    "            self.clf_ = {}  # 分组情况\n",
    "            for j in range(self.k_):\n",
    "                self.clf_[j] = []  # 每次迭代清空分组结果\n",
    "            for feature in data:\n",
    "                distances = []\n",
    "                for center in self.centers_:\n",
    "                    distances.append(np.linalg.norm(feature - self.centers_[center]))  # 欧式距离\n",
    "                classification = distances.index(min(distances))  # 单个数据的分组结果\n",
    "                self.clf_[classification].append(feature)  # 将单个数据添加到不同组中\n",
    "\n",
    "            print(f\"第{i+1}次迭代\")\n",
    "            print(\"中心点:\",end=' ')\n",
    "            pprint(self.centers_)\n",
    "            print(\"分组情况:\",)\n",
    "            pprint(self.clf_, width=120, indent=4, compact=True)\n",
    "            print(\"----------------------------------------------------------------------\")\n",
    "\n",
    "            prev_centers = dict(self.centers_)\n",
    "            for c in self.clf_:\n",
    "                self.centers_[c] = np.average(self.clf_[c], axis=0)  # 重新计算中心点坐标\n",
    "\n",
    "            # 中心点是否在误差范围\n",
    "            optimized = True\n",
    "            for center in self.centers_:\n",
    "                org_centers = prev_centers[center]  # 上一次的中心点坐标\n",
    "                cur_centers = self.centers_[center]  # 这一次的中心点坐标\n",
    "                if np.sum((cur_centers - org_centers) / (org_centers * 100.0 + 1e-6)) > self.tolerance_:\n",
    "                    optimized = False\n",
    "            if optimized:\n",
    "                break  # 两次中心点坐标比较相差无几后，结束循环\n",
    "\n",
    "    def predict(self, p_data):\n",
    "        \"\"\"  k均值预测数据  \"\"\"\n",
    "        distances = [np.linalg.norm(p_data - self.centers_[center]) for center in self.centers_]  # 欧式距离\n",
    "        index = distances.index(min(distances))   # 单个数据的分组结果\n",
    "        return index\n",
    "\n",
    "x = np.array([[0, 0], [1, 0], [0, 1], [1, 1], [2, 1], [1, 2], [2, 2], [3, 2], [6, 6], [7, 6],\n",
    "              [8, 6], [6, 7], [7, 7], [8, 7], [9, 7], [7, 8], [8, 8], [9, 8], [8, 9], [9, 9]])\n",
    "k_means = K_Means(k=2)  # 分成2类\n",
    "k_means.fit(x)\n",
    "for center in k_means.centers_:\n",
    "    pyplot.scatter(k_means.centers_[center][0], k_means.centers_[center][1], marker='*', s=150)  # 画出中心点\n",
    "\n",
    "for cat in k_means.clf_:\n",
    "    for point in k_means.clf_[cat]:\n",
    "        pyplot.scatter(point[0], point[1], c=('r' if cat == 0 else 'b'))  # 画出样本数据\n",
    "\n",
    "predict = [[2, 3], [6, 9]]\n",
    "for feature in predict:\n",
    "    cat = k_means.predict(feature)\n",
    "    pyplot.scatter(feature[0], feature[1], c=('r' if cat == 0 else 'b'), marker='x')  # 画出预测数据\n",
    "\n",
    "pyplot.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
